n_distinct(result$MA.id) # number of MAs
## [1] 23492
nrow(result) # 9 model result per each MA
## [1] 211428
table((result %>% distinct(MA.id, .keep_all = TRUE))$event)
## 
##     common infrequent       rare  very rare 
##      22074       1154        101        163
proportions(table((result %>% distinct(MA.id, .keep_all = TRUE))$event)) * 100
## 
##     common infrequent       rare  very rare 
## 93.9639026  4.9123106  0.4299336  0.6938532

Estimation

Analysis of heterogeneity estimate

There are three types of model that have \(\tau\)

  1. \(\tau\) from ivre, htevaguelogit and htelogit_tau50 (log ratios)

  1. \(\tau_a\) and \(\tau_c\) from ablogit and ablogit_tau50 (log odds)

  1. \(\tau_a\) and \(\tau_c\) from htebeta (prob)

Correlation of tau

cor(
  result_bayesian[result_bayesian$modelname == "htelogit_tau50",]$tau,
  result_bayesian[result_bayesian$modelname == "ablogit_tau50",]$taua
  )
## [1] 0.2477361
cor(
  result_bayesian[result_bayesian$modelname == "htelogit_tau50",]$tau,
  result_bayesian[result_bayesian$modelname == "ablogit_tau50",]$tauc
  )
## [1] 0.3046824

cor(
  result_bayesian[result_bayesian$modelname == "htebeta",]$taua,
  result_bayesian[result_bayesian$modelname == "ablogit_tau50",]$taua
  )
## [1] -0.08585503
cor(
  result_bayesian[result_bayesian$modelname == "htebeta",]$tauc,
  result_bayesian[result_bayesian$modelname == "ablogit_tau50",]$tauc
  )
## [1] -0.1115232

Analysis Significance Inference by model

modelname Sig_sum percentage_total
ivfe 10403 44.283160
ivre 9132 38.872808
ctevaguelogit 11553 49.178444
ctebeta 11402 48.535672
htevaguelogit 7308 31.108462
htelogit_tau50 6047 25.740678
htebeta 4363 18.572280
ablogit 8003 34.066916
ablogit_tau50 597 2.541291

table(result_sigsum$Sig_sum)
## 
##     0     1     2     3     4     5     6     7     8     9 
## 10947  1075   706  1407  1587   922   915  2312  3121   500
round(proportions(table(result_sigsum$Sig_sum))*100,3)
## 
##      0      1      2      3      4      5      6      7      8      9 
## 46.599  4.576  3.005  5.989  6.755  3.925  3.895  9.842 13.285  2.128
table(result_sigsum$event, result_sigsum$Sig_sum)
##             
##                 0    1    2    3    4    5    6    7    8    9
##   common     9948  977  632 1340 1524  902  883 2261 3109  498
##   infrequent  799   74   58   59   61   18   29   42   12    2
##   rare         65    8    7    5    2    2    3    9    0    0
##   very rare   135   16    9    3    0    0    0    0    0    0
round(proportions(table(result_sigsum$event, result_sigsum$Sig_sum), margin = 1)*100,2)
##             
##                  0     1     2     3     4     5     6     7     8     9
##   common     45.07  4.43  2.86  6.07  6.90  4.09  4.00 10.24 14.08  2.26
##   infrequent 69.24  6.41  5.03  5.11  5.29  1.56  2.51  3.64  1.04  0.17
##   rare       64.36  7.92  6.93  4.95  1.98  1.98  2.97  8.91  0.00  0.00
##   very rare  82.82  9.82  5.52  1.84  0.00  0.00  0.00  0.00  0.00  0.00

waic

For best model comparison among seven Bayesian models (exclude ivfe and ivre)

htelogit

result_htevaguelogit_tau = result %>% filter(modelname %in% c("ctevaguelogit", "htevaguelogit")) %>% group_by(MA.id) %>% mutate(diff_waic = diff(waic)) %>%  filter(modelname == "htevaguelogit") %>% ungroup()

mean(result_htevaguelogit_tau$diff_waic < 0) # htevaguelogit has better waic than ctevaguelogit
## [1] 0.9258471
result_htelogit_tau50_tau = result %>% filter(modelname %in% c("ctevaguelogit", "htelogit_tau50")) %>% group_by(MA.id) %>% mutate(diff_waic = diff(waic)) %>%  filter(modelname == "htelogit_tau50") %>% ungroup()

mean(result_htelogit_tau50_tau$diff_waic < 0) # htelogit_tau50 has better waic than ctevaguelogit
## [1] 0.919164
result_htelogit = result %>% filter(modelname %in% c("htevaguelogit", "htelogit_tau50")) %>% group_by(MA.id) %>% mutate(diff_waic = diff(waic)) %>%  filter(modelname == "htelogit_tau50") %>% ungroup()

mean(result_htelogit$diff_waic < 0) # htelogit_tau50 has better waic than htevaguelogit
## [1] 0.5440576

ablogit

result_ablogit_tau = result_narm %>% filter(modelname %in% c("ctevaguelogit", "ablogit")) %>% group_by(MA.id) %>% 
  mutate(diff_waic = diff(waic)) %>%  filter(modelname == "ablogit") %>% ungroup()

mean(result_ablogit_tau$diff_waic < 0) # ablogit has better waic than ctevaguelogit
## [1] 0.5141309
result_ablogit_tau50_tau = result_narm %>% 
  filter(modelname %in% c("ctevaguelogit", "ablogit_tau50")) %>% group_by(MA.id) %>% 
  mutate(diff_waic = diff(waic)) %>%  filter(modelname == "ablogit_tau50") %>% ungroup()

mean(result_ablogit_tau50_tau$diff_waic < 0) # ablogit_tau50 has better waic than ctevaguelogit
## [1] 0.6156923
result_ablogit_tau = result_narm %>% filter(modelname %in% c("ctebeta", "ablogit")) %>% group_by(MA.id) %>% 
  mutate(diff_waic = diff(waic)) %>%  filter(modelname == "ablogit") %>% ungroup()

mean(result_ablogit_tau$diff_waic < 0) # ablogit has better waic than ctebeta
## [1] 0.7459887
result_ablogit_tau50_tau = result_narm %>% 
  filter(modelname %in% c("ctebeta", "ablogit_tau50")) %>% group_by(MA.id) %>% 
  mutate(diff_waic = diff(waic)) %>%  filter(modelname == "ablogit_tau50") %>% ungroup()

mean(result_ablogit_tau50_tau$diff_waic < 0) # ablogit_tau50 has better waic than ctebeta
## [1] 0.7681851
result_ablogit = result_narm %>% 
  filter(modelname %in% c("ablogit", "ablogit_tau50")) %>% group_by(MA.id) %>% 
  mutate(diff_waic = diff(waic)) %>%  filter(modelname == "ablogit_tau50") %>% ungroup()

mean(result_ablogit$diff_waic < 0) # ablogit_tau50 has better waic than ablogit
## [1] 0.6549662
##                ctevaguelogit ctebeta htevaguelogit ablogit htebeta
## ctevaguelogit          0.000   0.741         0.074   0.486   0.280
## ctebeta                0.259   0.000         0.062   0.254   0.040
## htevaguelogit          0.926   0.938         0.000   0.920   0.754
## ablogit                0.514   0.746         0.080   0.000   0.265
## htebeta                0.720   0.960         0.246   0.735   0.000
## htelogit_tau50         0.919   0.945         0.541   0.925   0.782
## ablogit_tau50          0.616   0.768         0.199   0.655   0.443
##                htelogit_tau50 ablogit_tau50
## ctevaguelogit           0.081         0.384
## ctebeta                 0.055         0.232
## htevaguelogit           0.459         0.801
## ablogit                 0.075         0.345
## htebeta                 0.218         0.557
## htelogit_tau50          0.000         0.828
## ablogit_tau50           0.172         0.000
##                ctevaguelogit ctebeta htevaguelogit ablogit htebeta
## ctevaguelogit          0.000   0.000         0.000   0.000   0.000
## ctebeta                0.259   0.000         0.000   0.000   0.000
## htevaguelogit          0.926   0.938         0.000   0.000   0.000
## ablogit                0.514   0.746         0.080   0.000   0.000
## htebeta                0.720   0.960         0.246   0.735   0.000
## htelogit_tau50         0.919   0.945         0.541   0.925   0.782
## ablogit_tau50          0.616   0.768         0.199   0.655   0.443
##                htelogit_tau50 ablogit_tau50
## ctevaguelogit           0.000             0
## ctebeta                 0.000             0
## htevaguelogit           0.000             0
## ablogit                 0.000             0
## htebeta                 0.000             0
## htelogit_tau50          0.000             0
## ablogit_tau50           0.172             0

##                ctevaguelogit ctebeta htevaguelogit ablogit htebeta
## ctevaguelogit          0.000   0.133         0.006   0.584   0.075
## ctebeta                0.867   0.000         0.217   0.895   0.176
## htevaguelogit          0.994   0.783         0.000   0.958   0.559
## ablogit                0.416   0.105         0.042   0.000   0.062
## htebeta                0.925   0.824         0.441   0.938   0.000
## htelogit_tau50         0.991   0.832         0.551   0.961   0.672
## ablogit_tau50          0.575   0.180         0.070   0.650   0.114
##                htelogit_tau50 ablogit_tau50
## ctevaguelogit           0.009         0.425
## ctebeta                 0.168         0.820
## htevaguelogit           0.449         0.930
## ablogit                 0.039         0.350
## htebeta                 0.328         0.886
## htelogit_tau50          0.000         0.963
## ablogit_tau50           0.037         0.000
##                ctevaguelogit ctebeta htevaguelogit ablogit htebeta
## ctevaguelogit          0.000   0.000         0.000   0.000   0.000
## ctebeta                0.867   0.000         0.000   0.000   0.000
## htevaguelogit          0.994   0.783         0.000   0.000   0.000
## ablogit                0.416   0.105         0.042   0.000   0.000
## htebeta                0.925   0.824         0.441   0.938   0.000
## htelogit_tau50         0.991   0.832         0.551   0.961   0.672
## ablogit_tau50          0.575   0.180         0.070   0.650   0.114
##                htelogit_tau50 ablogit_tau50
## ctevaguelogit           0.000             0
## ctebeta                 0.000             0
## htevaguelogit           0.000             0
## ablogit                 0.000             0
## htebeta                 0.000             0
## htelogit_tau50          0.000             0
## ablogit_tau50           0.037             0